IS415 Web Blog: in-class-exercise-4

sf is mainly used for importing geospatial data
spatstat is used to perform spatial point patern analysis, might be a little outdated
maptools is used for map data conversion
ppp is the format supported by spatstat, hence we need maptools to convert spatial data from one type to ppp format
raster is used for handling raster-based geospatial data
tmap for mapping
tidyverse for non-spatial data

packages = c("maptools","sf","raster","spatstat","tmap","tidyverse")
 for (p in packages){
    if(!require(p, character.only = T)){
      install.packages(p)
  }
  library(p,character.only = T)
}

sg_sf = st_read(dsn = "data/shapefile",
                layer = "CostalOutline")

Reading layer `CostalOutline' from data source 
  `C:\JunLonggggg\junlong-is415\_posts\2021-09-06-in-class-exercise-4\data\shapefile' 
  using driver `ESRI Shapefile'
Simple feature collection with 60 features and 4 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 2663.926 ymin: 16357.98 xmax: 56047.79 ymax: 50244.03
Projected CRS: SVY21

mpsz_sf = st_read(dsn = "data/shapefile",
                layer ="MP14_SUBZONE_WEB_PL")

Reading layer `MP14_SUBZONE_WEB_PL' from data source 
  `C:\JunLonggggg\junlong-is415\_posts\2021-09-06-in-class-exercise-4\data\shapefile' 
  using driver `ESRI Shapefile'
Simple feature collection with 323 features and 15 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: 2667.538 ymin: 15748.72 xmax: 56396.44 ymax: 50256.33
Projected CRS: SVY21

This is a new read format,using rds file as input file type.

Childcare layer is just attribute table, aspatial data, no geometric data. Likewise for chas data.

We need to check the data even if the dataset come from the same agency, structure of a dataset may vary.

childcare = read_rds("data/rds/childcare.rds")
chas = read_rds("data/rds/CHAS.rds")

Converting from aspatial to geospatial

chas_sf = st_as_sf(chas, coords = c("X_COORDINATE","Y_COORDINATE"),
                   crs = 3414)

4326 is the EPSG for WGS84, so we need to set it first then convert the coordinate system.

SVY21 is in metres while WGS84 is in degrees

data scientist at data.gov usually handle geospatial data using kml or geojson file type

thus the API is introduced for a better structured geospatial data format

As a data analyst, we should choose the data that can benefit us the most when performing analysis

childcare$Lat = as.numeric(childcare$Lat)
childcare$Lng = as.numeric(childcare$Lng)
childcare_sf = st_as_sf(childcare, 
                        coords = c("Lng",
                                   "Lat"),
                        crs = 4326) %>% 
  st_transform(crs = 3414)

converting from sf format to Spatial* classes format

childcare = as_Spatial(childcare_sf)
chas = as_Spatial(chas_sf)
mpsz = as_Spatial(mpsz_sf)
sg = as_Spatial(sg_sf)

converting from Spatial* classes to sp format

SP layer do not have attribute/associate table

childcare_sp = as(childcare, "SpatialPoints")
chas_sp = as(chas, "SpatialPoints")
sg_sp = as(sg, "SpatialPolygons")

converting from sp to spatstat ppp format

Provides owin, the window study area, would be rectangle based on the xlim and ylim of the data

once converted, projection would be lost, when converting back to raster layer, we need to add in the projection

childcare_ppp = as(childcare_sp,"ppp")
chas_ppp = as(chas_sp,"ppp")

Keep interactive maps below 10 maps, we need to be aware of the limits web servers have

tmap_mode("view")
tm_shape(childcare_sf)+
  tm_dots(alpha = 0.4,
          col = "blue",
          size = 0.05) + 
tm_shape(chas_sf) + 
  tm_dots(alpha = 0.4,
          col = "red",
          size = 0.05)

converting tmap_mode back to plot

tmap_mode("plot")